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Abstract 

There exist a large literature on the application of (/-statistics to the out-of-equilibrium non- 
ergodic systems in which some degree of strong correlations exists. Here we study the distribution 
of first return times to zero, Pr( 0, t), of a random walk on the set of integers {0,1, 2,..., L} with a 
position dependent transition probability given by \n/L\ a . We find that for all values of a S [0,2] 
Pr(0, t) can be fitted by (/-exponentials, but only for a = 1 is Pji(0,t ) given exactly by a q- 
exponential in the limit L —>• oo. This is a remarkable result since the exact analytical solution of 
the corresponding continuum model represents Pr(0, t ) as a sum of Bessel functions with a smooth 
dependence on a from which we are unable to identify a = 1 as of special significance. However, 
from the high precision numerical iteration of the discrete Master Equation, we do verify that 
only for a = 1 is Pr(0, t) exactly a (/-exponential and that a tiny departure from this parameter 
value makes the distribution deviate from (/-exponential. Further research is certainly required to 
identify the reason for this result and also the applicability of (/-statistics and its domain. 
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I. INTRODUCTION 


It is well know that Boltzmann-Gibbs statistics is related to the exponential and Gaussian 
functions because these functions are the distributions corresponding to maximising the 
Boltzmann-Gibbs entropy. On the other hand, in cases where ergodicity is broken say 
due to strong correlations or long range interaction exist Boltzmann-Gibbs statistics can 
break down. It has been suggested that for a very wide class of such situations including 
experimental 1-3|, observational and examples of model systems ft] [l2 the so-called 
g-statistics is the appropriate generalization of Boltzmann-Gibbs statist ics fldl. Q . In this 
case ordinary exponential and Gaussian functions are replaced by g-exponentials and q- 
Gaussians, defined respectively as, 


exp u = 


[1 + (1 - q)u] l/{1 ~ q) , 1 + (1 - q)u > 0 


else 


( 1 ) 


and 


P(u) oc exp (— Bu 2 


( 2 ) 


which maximise Tsallis entropy S q = k (1 — JAp?) / (q — 1) under appropriate conditions 
jl3, 14], Much of the evidence for this suggestion comes from often very good quality fitting 
of the g-functions to simulated or observed data combined with the appeal of the idea that the 
q -functions are the maximum distributions an appropriate entropy of presumably very broad 
applicability. This view point suggests that when dealing with systems that takes us beyond 
the realm of Boltzmann-Gibbs, the g-functions should be the generic distribution functions. 
Accordingly we would expect g-functions describe the system whenever the choice of systems 
parameters leads to non-Boltzmann-Gibbs behavior. Given this, one would furthermore 
expect that even if we are not a priori able to determine the specific value of the g parameter 
in Eqns. dU) and ([2]) we should be able to obtain the relevant value of g by fitting since if q- 
statistics is the generic generalisation of Boltzmann-Gibbs statistics the family of g-functions 
should exhaust the functional forms distributions will assume. In order to understand better 
to which extend this attractive picture apply in reality, two of the authors together with C. 
Tsallis investigated a little while ago a one-dimensional Restricted Random Walk (RRW) 


model 


15| . The model consists of a random walker with a reflecting boundary condition 
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on the set of integers {0,1, 2, L}, where the probability to make a move is given by 
g{n ) = ( n/L ) a . The Master Equation for the model is straight forward to derive and to 
solve numerically by iteration to any desired numerical precision. We studied the sum of the 
positions of the walk after T steps. From very high precision numerics we concluded that in 
the limit L —> oo and for a certain scaling of T and L the distribution of this sum does not 
converge to the usual Gaussian given by the ordinary central limit theorem, but rather the 
limit distribution is a g-Gaussian. Surprisingly the g-Gaussian is only obtained for a = 1. 
For a 1 neither a Gaussian nor a g-Gaussian appears to describe the limit distribution of 
the sum of positions. It is not immediately clear why the value a = 1 for the exponent of 
the transition probabilities is a special case. 

In order investigate how a = 1 can be a singular value, we here investigate the distribution 
of first return times of the RRW. We solve the Master Equation for the probability to find 
the walker at position n at time step t analytically exactly in the continuum approximation 
and also solve numerically exactly the discrete version of the walker. From the analytic 
expression, which gives the distribution in terms of a Fourier-Bessel series, we are unable to 
identify anything singular about the value a = 1. However, the numerical solution indicates 
that indeed for a = 1 the return distribution is given by a g-exponential in the limit L —> oo, 
whereas for 1 the distribution is close but not exactly equal to any g-exponential. 

In the remainder of the paper we present firstly the analytic solution to the Master Equa¬ 
tion for the RRW model and investigate the Erst return distribution. We use the analytic 
solution presented here to investigate the nature of the dependence of the distributions on 
a, and we then assess the relevance of the g-statistics by numerically approaching the same 
problem. We close with a discussion of our results and their implications for how g-functions 
may relate to non-Boltzmann-Gibbs statistics. 


II. THE RESTRICTED RANDOM WALKER 

We consider a one-dimensional RRW model. We let n be the position of this model on 
the set {0,1, 2,..., L}, and s be the discrete time variable. The walker is confined to the 
integers between 0 and L, and therefore n G {0,1, 2..., L}. The motion of the walker is 
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controlled by the following process: 


where 


^s+1 


' ' 
n s + 1 with probability g(n)/2 

n s — 1 with probability g(n )/2 

n s with probability 1 — g{n) 

< " > 


( 3 ) 


9{n) 


a e (0, 2) 


( 4 ) 


Noting that the process for a normal random walker requires that g(n) = 1 so that the 
walker moves to the left or right with a probability of 1/2 at each time-step. 


III. THE DISCRETE MASTER EQUATION FOR THE DISTRIBUTION 

The Master Equation for the distribution is given by 
P n (n, s + 1) =P n (n, s ) + ^g(n - 1) 

x P n {n - l,s) T \g( n + 1) 
x P n (n + 1, s) - g(n)P n (n , s). 

We note that the Master Equation given in Eq. (J5]) is valid for the bulk sites. At n = L 
we have a reflective boundary condition, while at n = 0 we have an absorbing boundary 
condition. The probability of the return for the RRW at n — 0, is therefore as follows: 

Pr d (^i s + 1) = P n ( 0, s + 1) — P n (0, s) (6) 


IV. THE CONTINUOUS EQUATION FOR THE DISTRIBUTION 
A. The Master Equation in continuous form 


In order to obtain the continuous form of the Master Equation we introduce Ax = —, 

r 12 1 

i.e., x = nAx G < 0,1 j>. Similarly we scale time and choose a time increment 
1 s 

At = — or t = —. This ensures that the continuum approximation becomes exact in the 

1j Ju 

limit L —> oo. Finally we replace g{n) by g(x), where g(x) = \x\ a . 
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As L —» oo, the Master Equation takes the form of the following continuum diffusion 
equation 

rtPlr A 1 r ) 2 

(7) 


dP (*’*) = J_ u- t) i 


A similar problem was studied in 16J. Using the separation of variables the solution for 
P(x,t ) can be found as 


P(x,t) = ^2 A x P x (x) exp tj ■ 

Assuming Q\{x) = x a P\(x), the differential equation in space becomes 


A 2 


( 8 ) 


d tQx( x ) + = o. 


(9) 


dx 2 A v ' Ax a 

The general solution of this equation can be written in terms of Bessel and Neumann func¬ 
tions |d, p. 362, 9.1.51]: 


Q\( x ) =A yfx J\/( 2 -a) 

+ B y/x Y\ j ( 2 -a) 


A x l ~ a B 
2 — a 
Ax 1 -"/ 2 
2 — a 


( 10 ) 


The boundary conditions are the following: At x = 0 given that P\{x) is finite, Q\(x = 0) 

is zero. At x = 1 we have a reflective boundary, therefore the current J, at x = 1 is zero i.e. 

dP dJ 

J(x,t) L=i= 0. Further we know that —— = ——. Using this fact together with Eq. (17]) 

at ox 

\ c) d 

we can conclude that J(x,t) = - - 7 — \x a P(x,t) ]. We know that if — Q\(x) L=i= 0, this 

2 L ox dx 

condition is satisfied. 

The boundary condition at x = 0 implies that B — 0. The boundary condition at x = 1 


implies that 


A 


2 — a 


has to be the zero root of the Bessel function of order 


Ar 


2 — a 


Therefore - = j\/( 2 - a ) -1 n, where j 1 /( 2 - a ) -1 n is the ^th zero root of the Bessel function 

2 — a ’ 

1 


of order 


2 — a 


- 1 


So we arrive at the follwoing solution for P(x, t ) 

J \/ ( 2 - a) (j 1 /( 2 -a) -1 ,n x 


B( X A) = ^Z A r. 


l-a/ 2 ^ 


n= 1 


x a- 1/2 
\ 2 .'2 


x exp 


(2 CL) J 1/(2-a) -1 
8 L 


( 11 ) 
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For the special case of a = 1, it turns out that 


P(x,t) = ^A n 

n=l 



exp 


do,n A 
8 L / 


( 12 ) 


B. The solution of the coefficient A n 

Our initial condition for the problem is that at t — 0, all walkers are placed at n — 1. 
We attempt to solve coefficient A n for this initial condition. At t — 0 the distribution is as 
follows: 


P{x, 0)=5> r 


Jl/(2-a) (jl/(2-a) -1 ,n X 


l-a/2j 


n= 1 


X 


a- 1/2 


(13) 


To solve coefficients A„ we refer to the Fourier-Bessel series 


18 


, Chap. 8]. We note that 


we have a slightly different problem here. We therefore present a modified version of the 
Fourier-Bessel series below. 


We note that Eq. (J9J) is a standard Sturn-Liouville equation pj9|, and as a result its 
eigenfunctions are orthogonal. Therefore if a and b are different zeros of the Bessel function 
of order u — 1, then f J u (ax) J u {bx ) xdx = 0. If a —> b the solution to this integral is J%(a )/2 
[20I ] . Hence one can easily write 


Jv{jv—l,m,X)Ji/{jv—l,nX')xdx \^u (.A—l,m)] 


where S m ^ n is the Dirac delta function. 


(14) 


Using the above result, we obtain 

2(1-a/2) 


A n — 




1 (il/(2 -a) -1 ,n) "s/r 

For the special case of a = 1, this result reads 


Ji 


1 /( 2 -“) 


Jl/(2-a) —1 ,n 

L l ~ a / 2 


(15) 


An — 


:Jl\ 


J0,r 


(16) 


) vV \v / I / 

We now consider the asymptotic behaviour for large t. Using the asymptotic expansion 
for J u (z) 3, P- 199, 7.21], as n —> 00 , A n is asymptotically equal to 
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Aj ~'\/^ 7r ^1/(2 -a) —l,n (l — g) ^ a ^ 

Jl/(2-a) — 1 ,n 7T 

x cos 1 L i-a/2 2 (2-a) 



( 17 ) 


This implies that as long as — is not too small, the higher terms in the sum in Eq. ([IT]) can 

1j 

be ignored and P(x, t ) will decay exponentially for t > L or in terms of the discrete time 
step variable s > L 2 . 


V. ANALYTIC SOLUTION FOR THE PROBABILITY OF RETURN 


Since the boundary at x - 
1 d 

current J[x,t) = — — — \x 
2L ox 


= 0 is absorbing the first return time probability is given by the 
“ P(x,t )] evaluated at x — 1/L and we find 


\[L 00 

Pft(0> ^ ^ A n J\/(2-a ) ( Jl/(2-a) -1 ,r 


n =1 


x exp 


(2 a) J 1/(2-a) -1 
8L 


^ l-a/2 


t . 


( 18 ) 


Fig. m show a comparison between the exact numerical solution to the Master Equation 
in Eq. (J5]) and the analytic solution in Eq. f[T8]) . The number of terms is fixed to N = 1000 
which ensures that the plotted analytic result is correct within an accuracy smaller than the 
width of the plotted graph. 

Unfortunately we have not been able to reduce the series in [18] to a simple compact 
expression. In particular we were unable to relate analytically .Pr( 0, i) to the ^-exponential 
forms given in Eq. ([T]). Since the Bessel functions in the series for Pr( 0,£) are smooth 
functions of the index 1/(2 — a) and 1/(2 — a) — 1 one will not expect that the description in 
terms of (/-exponentials would be more appropriate for certain values of the exponent a in 
Eq. (J3]) than for others. However, as well see in the next section high precision numerical 
analysis suggests otherwise. 


VI. NUMERICAL ANALYSIS USING (/-EXPONENTIALS 

Next we solve by numerical iteration the discrete Master equation Eq. (J5]) to obtain 
numerically exact results for the return time distribution in Eq. (jUJ) . These results are exact 
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L 

normalization 

20000 

0.9899 

50000 

0.9899 

100000 

0.9899 

200000 

0.9900 


TABLE I. For various values of L with a = 1, the normalizations attained from simulations are 
given. 


within the numerical precision of the computer. 

Using the definition of (/-exponentials in Eq. (IT]) , one can easily write the normalized first 
return distributions as 


Pr d (0, s) = P Rd (0, 0) exp ? (~P q s) + (19) 

where j3 q > 0 is the only free parameter to be fitted since the value of q can be obtained from 
the asymptotic behaviour of the slope, which must be equal to l/{q— 1). The comparison 
between the discrete analytic results (Eq. ((5])) and (/-exponential (Eq. ([19)) ) is given in Fig. [2] 

Because the numerical study is limited to finite time the normalization will not be exact, 
since we lack some contribution to the return distribution from the tail of very large times. 
We therefore check (see Table 1) the extend to which the normalization is fulfilled for the 
return times we are able to handle. 

A very good agreement can easily be seen from the figure, but at this point we make one 
more step to quantify this agreement better. To achieve this, we define a quantity A, which 
is the area between the curve of the exact discrete result and related (/-exponential function. 
We know that for finite L the return time distribution cuts off exponential for (discrete) 
time larger than L 2 . This means that the (/-exponential form will only be valid for all times 
in the limit L —> oo. If the distribution, in this limit, is exactly a (/-exponential, then A 
must approach zero as the system size tends to infinity. 

We plot in Fig. |3] this quantity as a function of system size for a = 1 case and also 
for two very tiny departures from this case. We note that for a = 1 the limit L oo is 
indeed consistent with a (/-exponential. For a ^ 1 by just very small amount we obtain 
lini/^oo A > 0 and therefore strictly speaking the return distribution is not equal to a 










g-exponential for a ^ 1 . 


VII. DISCUSSION AND CONCLUSION 


We have analytically and numerically analysed a simple random walk with a position 
dependent transition rate with the aim to better understand when and why (/-functional 
forms describe the statistical properties of non-conventional systems. We focus on the first 
return time distribution. 

The analytic solution suggests that the return distributions are smoothly dependent on 
the a parameter, which determines the transition probabilities of the RRW model. But de¬ 
spite of being unable analytically to relate the return time distribution to (/-exponentials, we 
found numerically that, in the limit of large systems sizes, the (/-exponential fits remarkably 
well the numerical exact result. And indeed for a special value a = 1 the numerical analysis 
demonstrate that the return distribution becomes equal to the (/-exponential with q = 3/2 
when the system size is taken to infinity. 

The question then is what is so special about a = 1 that for this value of a the q- 
exponential is equal to the considered distribution and not just a very good fit? 

One can identify a few ways in which a = 1 is special. The mean of the return time can 
be calculated directly (see e.g. 21 ] and 22 ]) without the need of knowing the entire return 
time distribution. The average of the number of discrete time steps for the first return from 
position n = 1 to n — 0, which we denote by (sr), is to leading order in the system size L 
given by 


L/(l — a) for 0 < a < 1 
( Sr) — \ L In L for a = 1 

L a / (a — 2) for 1 < a < 2. 

This indicates one way in which a = 1 is marginal. 

The (/-exponential found for a = 1 has index q — 3/2 which corresponds to a functional 
form f(x) = a/( 1+x) 2 . We want to mention that this form is also the functional dependence 
of the distribution of extinction times for a critical birth-death process in which the rate of 
death is equal to the rate of birth [23j. This result is likely to be related to the RRW studied 
here when a — 1, since the birth-death process is related to a critical branching process, 
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which on the other hand can 
proportional to the position 


ae related to a random walk with a transition probability 


24] 


Given these remarks it is difficult to decide whether the success of the fit to g-exponentials 
is caused by some fundamental underlying principle or is accidental in origin. 
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FIG. 1. The comparison between the discrete Master Equation (Eq. ([5]) ) and the continuous 
analytic solution (Eq. (1181) 1 for the first retrun time (as function of discrete time)for (a) a = 0.75, 
(b) a = 1 and (c) a = 1.25. 
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FIG. 2. The comparison between the discrete Master Equation (Eq. ([5]) ) and (/-exponential func¬ 
tions for the normalized first return distributions (as function of discrete time) for the cases (a) 
a = 0.75, (b) a = 1 and (c) a = 1.25. 
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FIG. 3. The quantity A as a function of the system size. 
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